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ABSTRACT 

The cosmic infrared background (CIB) radiation is the cosmic repository for energy 
release throughout the history of the universe. The spatial fluctuations of the CIB 
resulting from galaxy clustering are expected to be at least a few percent on scales 
of a degree, depending on the luminosity and clustering history of the early universe. 
Using the all-sky data from the COBE DIRBE instrument at wavelengths 1.25 - 100 
jum we attempt to measure the CIB fluctuations. In the near-IR, foreground emission 
is dominated by small scale structure due to stars in the Galaxy. There we find a 
strong correlation between the amplitude of the fluctuations and Galactic latitude after 
removing bright foreground stars. Using data outside the Galactic plane > 20°) and 
away from the center (90° < I < 270°) we extrapolate the amplitude of the fluctuations 
to cosec|6| = 0. We find a positive intercept of 6F,ras = 15.51?;^, 5.91^;^, 2.4l:gJ, 2.0lg:f 
nWm~^sr~^ at 1.25, 2.2, 3.5 and 4.9 /xm respectively, where the errors are the range of 
92% confidence limits. For color subtracted maps between band 1 and 2 we find the 
isotropic part of the fluctuations at T.Glgi nWm~^sr~^ . Based on detailed numerical 
and analytic models, this residual is not likely to originate from the Galaxy, our 
clipping algorithm, or instrumental noise. We demonstrate that the residuals from the 
fit used in the extrapolation are distributed isotropically and suggest that this extra 
variance may result from structure in the CIB. We also obtain a positive intercept from 
a linear combination of maps at 1.25 and 2.2 fim. For 2° < 9 < 15°, a power-spectrum 
analysis yields firm upper limits of (9/h°) x 5F^^s(d) < 6, 2.5, 0.8, 0.5 nWm~^sr~^ at 
1.25, 2.2, 3.5 and 4.9 fim respectively. From 10 - 100 fim, the dominant foregrounds 
are emission by dust in the Solar system and the Galaxy. There the upper limits on 
the CIB fiuctuations are below 1 nWm~^sr~^ and are lowest (< 0.5 nWm~^sr~^ ) at 
25 fxm. 

Subject headings: Cosmology - Cosmic Background Radiation - Galaxies: Clustering - 
Galaxies: Evolution 
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1. Introduction 

Diffuse backgrounds and their structure contain important information about tfie fiistory of 
the early Universe, when discrete objects either do not exist or are not accessible to telescopic 
studies. The formation and early evolution of galaxies should have generated radiation redshifted 
into the infrared bands (e.g. Partridge and Peebles 1967; Bond et al. 1986). This cosmic infrared 
background (CIB) may come from the entire history of the Universe between the epoch of last 
scattering, mapped by the microwave background, and the present day. 

The CODE Diffuse InfraRed Background Experiment (DIRBE) (Boggcss ct al. 1992; 
Silverberg et al. 1993; Weiland et al., DIRBE Explanatory Supplement, 1998) mapped the entire 
sky with 0.3° pixels and 0.7° resolution in ten bands between 1.25 and 240 jum. The results have 
been published (Hauser et al. 1998; Kelsall et al. 1998; Arendt et al. 1998; Dwek et al., 1998), 
based on a zodiacal light model believed accurate within a few percent (Kelsall et al. 1998). Their 
Galactic dust emission model (Arendt et al. 1998) is derived from fits to the DIRBE data and 
comparison to hydrogen maps, and accounts for the variation of dust temperature. Their model 
for Galactic starlight is derived from external star count models, with no free parameters. The 
residuals from the modeling were not significantly above the uncertainties except at 140 and 
240 /Ltm. The FIRAS instrument on COBE gave similar answers, (Fixscn, Mather, Bennett, and 
Shafer, 1998), using three different methods for removing the foregrounds. The agreement is 
significant because the instruments were calibrated independently. The measured values of the 
far-IR background are comparable to the total observed Galactic emission at visible and near-IR 
wavelengths, and imply that about half of the luminosity of the universe has been obscured by 
dust and converted into far-IR radiation. Some distant galaxies might be more reddened or 
obscured by dust than otherwise expected, and the CIB and its fluctuations might be brighter 
than predicted from visible band galaxy counts. 

In this paper we continue our search for a CIB by trying to detect its spatial structure in 
the DIRBE data. The method is similar to that first suggested by Gunn and later applied in 
the visible (Shectman 1973, 1974) and UV (Martin and Bowyer 1989). Recently, Vogeley (1998) 
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applied the method to the Hubble Deep Field constraining the amount of the visible cosmic 
background, concluding that there is very little visible cosmic background from unresolved sources. 
They considered arcsecond scales where fluctuations are dominated by random parts of individual 
galaxy (Vogeley 1998) or cluster (Shectman 1974) profiles. We apply the method to degree scales 
where the dominant contribution comes from galaxy clustering. Our previous results are described 
by Kashlinsky, Mather, Odenwald and Hauser (1996a; Paper I) and Kashlinsky, Mather and 
Odenwald (1996b; Paper II). 

If we assume that the galaxies producing the GIB are clustered, the amplitude of the resulting 
GIB 2-point correlation function, C{9), depends on the rate of the flux emission, dF/dz, and 
on the galaxy 2-point correlation function, S^{r;z). If the latter is known on the linear scales 
subtended by a given angle, measurements or upper limits on C{9) can constrain the levels of 
the GIB emitted by clustered material. Measuring or limiting the structure of the GIB can be 
especially valuable in the mid- to far-IR bands where the foreground emission is very bright, but 
smooth. 

The plan of the paper is as follows: In Sec. Q we define the quantities used and provide 
mathematical background. We show that for the measured galaxy correlation function we expect 
fluctuations of about 5 to 10% in the GIB flux on 0.5° scales, almost independently of the 
particular mechanism of the GIB production. We discuss our search for the GIB structure in the 
DIRBE data between 1 and 100 ^m. The last two bands of DIRBE - 140 and 240 ^um - are not 
useful for this because of the large instrument noise there. Section |3| deals with the beam profile 
and the spatial window function for power spectrum analysis. Section ^ discusses the data for 
A < 5 fim, where the foregrounds are dominated by stellar emission from the Galaxy. We find 
a significant trend in the amplitude of the measured fluctuations with Galactic latitude. When 
our star removal algorithm is applied to simulated data, we recover the measured slope of the 
correlation with Galactic latitude, and extrapolation of the simulated fluctuations to cosec|6| = 
leads to a zero intercept. We derive a formula relating the measured fluctuations to star counts 
and show that for a plane-parallel Galaxy there should be a simple relationship between the b 
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dependence of the fluctuations and the slope of the star counts function at the Galactic pole. 
Then we show that for 90° < / < 270° and |6| > 20° the measured structure of the Galaxy flts the 
plane-parallel model. Furthermore, the scatter in the C(0) - cosec|6| plot is low enough to allow 
an extrapolation to cosec|6| = 0, with positive intercepts in all four near-IR bands. We perform 
a power spectrum analysis but most of the measured structure is due to fluctuating star counts. 
If the star removal algorithm removes stars fainter than 5 times the confusion noise limit, there 
are too few pixels left to calculate a power spectrum. Nevertheless, limits are sufficiently low to 
be of interest in testing theoretical predictions. We test our methods with numerical and analytic 
models of the Galaxy and possible GIB contributions and find good agreement with the data. 
Results from the mid- to far-IR bands, where foreground emission structure is dominated by the 
large scale gradients from dust emission in the Solar System and the Galaxy, are presented in Sec. 
^. The upper limits on these fiuctuations are < 1 nWm~^sr^^ , depending on wavelength, and 
provide much stronger tests of the predictions of galaxy counts than do direct measurements from 
inside the zodiacal dust cloud. We summarize our results in Sec. |6[ 

2. Theoretical preliminaries 

Here we provide a general mathematical basis without specific cosmological and galaxy 
evolution models. We will use the data to constrain the GIB properties in a model-independent 
way, and then to constrain models of galaxy evolution. We extend our previous work (Paper 1,1; 
Jimenez & Kashlinsky 1999) to show the power spectrum of the fiuctuations and to estimate the 
typical amplitude of the fiuctuations in a model-independent way. 

We start with definitions. The surface brightness in the GIB per unit wavelength will be 
denoted as Ix, per unit frequency as 1^, and per logarithmic wavelength interval F = A/a = i^Iu, 
and we call them all "flux." The fluctuation in the GIB flux is then 6F(x) = F{x) — (F), where 
X is the two dimensional coordinate on the sky and (...) denotes ensemble averaging. We also use 
the two-dimensional Fourier transform, 5F{9) = (27r)^^/ 5Fgexp(— iq • 6)(Pq. 
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If SF{x) field is a random variable, then it can be described by the moments of its probability 
distribution function. The first non-trivial moment is the projected 2-dimensional correlation 
function C{6) = {5F{x + 6)SF{x)). The 2-dimensional power spectrum is ^2(3) = (I'^-^gP)) where 
the average is performed over all phases. The correlation function and the power spectrum are a 
pair of 2-dimensional Fourier transforms and for isotropically distributed signal are related by 

C{.e) = ^j^ P2{q)Mqe)qdq, (1) 

P2{q) = 27r / c{e)Jo{qe)ede, (2) 

Jo 

where Jn{x) is the n-th order cylindrical Bessel function. If the phases are random, then the 
distribution of the flux field is Gaussian and the correlation function, or its Fourier transform the 
power spectrum, uniquely describe its statistics. In measurements with a finite beam radius 'Q 
the intrinsic power spectrum is multiplied by the window function W of the instrument. Another 
useful quantity is the mean square fluctuation within a finite beam, or zero-lag correlation signal, 
which is related to the power spectrum by 

{{SFfU = — / P2{q)WTH{q^)qdq 

~ 2^^^^2(g)|g~7r/2tf. (3) 

For a top-hat beam the window function is Wth = [2Ji(a:)/a:]^ = 0.5 at a; 2± 7r/2 where x = q-d, 
and hence the values of q~^ correspond to fluctuations on angular scales of diameter ~ 7r/g. 

The GIB flux and its structure are measured in projection on the celestial sphere and reflect 
both the 3-dimensional clustering pattern of the galaxy distribution and the rate of emission at 
redshift z. We introduce the 3-dimensional two-point correlation function of galaxy clustering, 
^(r), and its 3-dimensional power spectrum, Pz{k). These are related via 3-dimensional Fourier 
transforms, and assuming isotropy are related by 

1 

i{r) = P3{k)jo{kr)k^dk. (4) 

The projected GIB correlation function is related to the underlying two-point correlation 
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function of the galaxy distribution and the rate of the CIB flux emission via the Limber equation 

r r fIF dF 

^(^)= / / :^:^e(ri2;^)d^id^2, (5) 

Jzj_ Jz2 dzi dZ2 

where ri2 is the proper length subtended by the angle 6 and redshifts zi,Z2- In the limit of small 
angles, ^ ^ 1, the Limber equation becomes (e.g. Peebles 1980) 

C{d)= / dz( — ) / ^{nr,z)du. (6) 

Jo \dz J J-oo 

For Robertson- Walker metrics the proper separation is given by 

, , dx/dz x'^(z)9'^ 

= <7T=i?"' + (Trip 

= c'{^j-V + d\{z)e'. (8) 

Here dA{z) is the angular diameter distance, x{z) = dA{z)/{l + z), u is the integration variable 
and t is the cosmic time. 

The comoving volume occupied by a unit solid angle in the redshift interval dz is 
dV/dz = (1 + z)x'^{z)cdt/dz, and the power received per unit wavelength and collecting 
area in band A from each galaxy with absolute bolometric luminosity L at redshift z is 
[L/(47rx^(l + z)^)]fx{i^', z). Here f\d\ is the fraction of the total light emitted in the wavelength 
interval [A; A + dX\ and the extra factor of (1 + z) in the denominator accounts for the fact that 
the flux received in band A comes from a redshifted galaxy. The contribution to the total CIB flux 
from the redshift interval dz is given by 

dF Rh 1 d{Hot) ^ A 

= 47 (TT^^i^ ^ A(z)[AA,(^; z)], (9) 

where the sum is taken over all galaxy populations contributing flux in the observer rest-frame 
band at A, and fx characterizes the spectral energy distribution (SED) of galaxy population i. 
Here Rh = cHq^ and C{z) = J ^{L; z)LdL is the comoving luminosity density from galaxies with 
the luminosity function ^{L; z) at epoch z. 

It is illustrative to study the redshift dependence of the flux rate production, Eq. (^). At 
small redshifts the factor (1 -|- z)^'^dt/dz varies little with and the rate of flux production is 
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governed by the comoving bolometric luminosity density C{z) and the SED of the galaxy emission 
f\. If the luminosity evolution at these redshifts is small, the rate of flux production is governed by 
the SED shape. If /a(A) increases towards shorter wavelengths then dF/dz increases with z. For 
\f\ = const, and no luminosity evolution, the rate is roughly constant with small z. In addition, 
there is observational evidence for an increase of C{z) out to z ~ 1 in the galaxy rest-frame UV 
to near-IR (1 ^um) bands (Lilly et al. 1996). At sufficiently high redshifts, the evolution in Eq. (P) 
would be offset by the factor (1 + z)~'^dt/dz, so that the rate of production would be cut off at 
sufficiently large z. This factor is responsible for resolving Olbers' paradox even for a flat SED. 

The SED for rest-frame A < lO/^m is dominated by stellar emission, with a peak at visible 
wavelengths and a decrease for A > 0.7 iJ,m. Consequently, most of the predicted J band CIB 
comes from redshifts z ~ 0.3 — 1, which shifts the visible emission of normal stellar populations 
into the J band; cf. Yoshii and Takahara (1988). In the M band at 5 ^m, most of the predicted 
CIB comes from z > 1 — 2. At A > 10 /xm, the emission is dominated by galactic dust and the 
situation is reversed, so f\ increases with wavelength roughly as A" with a ~ 1.5. Hence, the dusty 
star-burst galaxies observed by IRAS at low redshifts should make the dominant contribution to 
the 10 fj,m CIB. In the far-IR, the measured CIB found by DIRBE can have large contributions 
from high redshifts. 

Measurements of the correlation function C{6) = {6F{x + 0)5F{x)) have an important 
advantage over direct determinations of the CIB spatial power spectrum, because they are 
immune to the discontinuities in the maps created by point source removal. (Because, in 
practice, C{6) is evaluated from masked data, the uncertainties of the C{9) points are themselves 
correlated.) However, interpreting C{6) in terms of dF/dz and ^(r) can be cumbersome because 
the right-hand-side of (P) contains a double integral and the 3-dimensional correlation function is 
not always positive. By definition, ^(r) must be negative on large scales so that /q°° ^(r)r^dr = 0. 

A simpler method is to work with the Fourier transform of the correlation function, P2{(l)^ 
which contains the same information as C{9) and is easier to interpret. As for C{6) the effects 
of the mask will produce correlations of the power spectrum points with each other. Also, the 
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measured power spectrum is the convolution of the Fourier transform of the mask with the true 
power. Although the uncertainties in C{9) and ^2(9) can be evaluated from theory, we evaluated 
them from the data by comparing multiple measurements. 

The Limber equation can be rewritten directly in terms of the power spectra, substituting 



Eq. (^) into (^) and using J^^jo{y/x'^ + y'^)dy = ■kJq{x) to obtain 

^ ^ 27rio \dz J cdt/dz 

P3{k;z)Jo{kdA9)kdk. (10) 

Substituting this expression into (^) and using the orthogonality relations for the Bessel functions, 
/o°° Jo{a9)Jo{/39)9d9 = a-^6D{a - (5), leads to 

Equation(^) involves only one integration and the kernel is always positive. Hence it may be 
easier to derive cosmologically interesting quantities such as F and ^(r) or P2,{k). Determining 
P2{q) directly from the data is difficult owing to the masking created by point source removal, and 
it seems best to compute C{9) and then derive P2{q) from it. 

Figure [l|a plots the argument, q{l + z)/x{z), of P3 in the right-hand-side of Eq. (|l^) for the 
largest wavenumbers (smallest scales) probed by DIRBE at various redshifts z. The dotted line 
corresponds to = 1, J7a = 0; the solid to 17 = 0.1, = and the dashes to O = 0.1, J7a = 0.9, 
where Q\ = 3Hq{1 — Q) denotes the contribution of the cosmological constant. For redshifts 
contributing most of the CIB flux (z >0. 1-0.2, cf. Paper I), the DIRBE instrument probes scales 
which are in the quasi-linear or linear regime, and which thus can be approximated as having 
evolved at the same rate with time. Thus for DIRBE scales, one can write Psik; z) ~ Psik; 0)^'^(2;), 
(e.g. Peebles 1980), where ^^(z) accounts for the evolution of the clustering and is normalized to 
^'2(0) = 1. Then Eq. (0) can be rewritten in a more compact way: 
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with 



H -= 1 n3^ 

{i + zy^i + nz + nA[{i + zy^-i\' ^ ' 



The left-hand-side of (12) is the same order of magnitude as the mean square fluctuation (|^, and 
we have defined 

A\k) = R-H^k^P-i{k;Q). (14) 

The quantity A(A;) is roughly the fluctuation over the line-of-sight cylinder of length Rh and 
diameter A:^^. For relevant scales and spectra of density fluctuations, A(A;) increases with k. 
Hence, the CIB fluctuation on a scale of /c"^ is ~ ^/WP^X^ ~ F^ik/Rn)- 

The present day power spectrum of galaxy clustering has been measured on scales 
corresponding to at least the smallest angular scales probed by DIRBE. The most accurate 
measurement comes from the APM survey data on the projected angular correlation function 
(Maddox et al. 1990). Baugh and Efstathiou (1993) deprojected the APM data to obtain the 
underlying power spectrum of galaxy clustering, P^{k). Kashlinsky (1998) used the current large 
scale data, and the abundance of objects at high redshifts, to reconstruct the pregalactic density 
field over six orders of magnitude in mass. His result requires significant fluctuation power on 
small scales, an early epoch of galaxy formation, and high levels of the CIB and its fluctuations. 
Fig. |l]b plots the data from the Baugh and Efstathiou (1993) deprojection. On small scales (large 



A;), A (A;) oc k , so the integrand in Eq. (12) behaves as oc z with an integrable singularity 



at z ^ 0. While a non-negligible part of the clustering part of the CIB fluctuations comes from 
nearby galaxies, much of it arises from galaxies at z > 1. 

In Papers I and II we estimated ~ 5-10% CIB fluctuations on a scale of 0.5°. One can also 
see this in a more intuitive way from Eq. (|l2|). Fig. |l|a plots q{l + z)/x{z), the largest wavenumber 
that enters on the right hand side of (12), for q^^ = 0.5°. In the Friedman-Robertson- Walker 



Universe, q/dA{z) reaches a minimum at z ~ 1 — 3, and its value at the minimum depends weakly 
on cosmological parameters. Thus ( |T2|) can be rewritten as an inequality: 

.(..,^,)/;(-)^g.. 
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On linear and quasi-linear scales over the range of redshifts that contribute most to the GIB, the 
quantity ^'^(z)/-^^ depends weakly on z (e.g. on linear scales = (1 + z)^^ if = 1). The 
integral is of the same order of magnitude as F'^ because the term dF/dz is a peaked function 
with a full width at half maximum of order unity. Hence, the relative fluctuation in the GIB is 
~ A(min[^|-^^^]). Fig. |^b plots A (A;) vs. k using the power spectrum of APM galaxies from Baugh 
& Efstathiou (1993); it is an increasing function of k. At the wavenumber of the minimum in Fig. 
||a, the GIB fluctuation can be as large as ~ (5-10)%. In models where the bulk of the GIB comes 
from higher redshifts, this number may be somewhat smaller. Because this method measures a 
two-point process, the results constrain a measure related to the mean square of the GIB emission 
rate, i.e. ~ / dz. 

Many models have calculated the expected GIB over the 1 - 100 fim range (e.g. Partridge 
and Peebles 1967; Stecker et al. 1977; Bond et al. 1986; Fah et al. 1996; Wang 1991; Beichman 
and Helou 1991; Franceschini et al. 1991; Gole et al. 1992; Malkan and Stecker, 1997; Jimenez 
and Kashlinsky 1999; Dwek et al. 1998). The models are normalized to galaxy counts and predict 
a typical flux of F ~ 10 nWm^^sr^^ . Deep K band counts of galaxies (e.g. Gowie et al. 1994; 
Djorgovski et al. 1995) or at 12 - 100 fim (Hacking and Soifer 1991) suggest minimal fluxes of at 
least a few nWm~^sr~^ . Therefore the GIB fluctuations on 0.5° scales may be ~ 1 nWm~^sr~^ . 

The APM measurements of the galaxy correlation function in the blue band might not 
apply to the infrared. However, the r band Palomar survey (Picard 1991), and the R-band Las 
Gampanas survey (Shectman et al. 1996, Lin et al. 1996) give results identical to the APM 
survey. On the other hand, for IRAS galaxies on small (non-linear) scales the correlation function 
has a lower coherence length. Saunders et al. (1992) show that on small scales IRAS galaxies 
at 60/xm cluster with S,iras = (r/r^jRAs)^^^^'^'^'^ where r^jRAS — 4/i~-'^Mpc and Jiras — 0.6, 
as opposed to r^: = 5.5/i~^Mpc and 7 ~ 0.7 for the APM galaxies (Moore et al. 1994). The 
slight difference could be due a tendency of the IRAS dusty star-burst galaxies to avoid the 
central regions of rich clusters of galaxies. However, on larger (linear) scales where the galaxy 
clustering pattern presumably traces the pregalactic density field, IRAS galaxies as measured by 
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the QDOT counts-in-cells analysis (Saunders et al. 1990) are consistent with the APM galaxies 
power spectrum. These linear to quasi-linear scales are relevant for the DIRBE beam size, so we 
use the APM numbers for all infrared wavelengths. In surveys with smaller beams, the differences 
will be more pronounced, but the differences between the IRAS and APM correlation functions do 
not lead to appreciably smaller mid- to far-IR CIB fluctuations. 

In addition to the clustering term, there is a white-noise component in the correlation signal 
due to individual galaxies (e.g. Peebles 1980). Its amplitude at zero lag is 

Rh f J^{L;z)L'^dLHodt 



Cwn(O) — 6 



^^^'"167r2y x2(z)(l + z)8 dz 

dz. (16) 



A ^2 



Like white noise from discrete stars, the white noise from galaxies is dominated by the nearest 
objects, unlike the clustered component (P), because of the presence of x'^{z) in the denominator of 
the integrand in (|l6|). Because the galaxies are almost undetectable against the star fluctuations, 
the galaxy white noise is also negligible. 



3. DIRBE data and beam profile 

The COBE DIRBE instrument provided an all-sky 41 week survey with a ten-band photometer 
(Boggess et al. 1992). The DIRBE bolometer measurements at 140 and 240 fim are too noisy for 
our purposes. The remaining eight bands are centered on wavelengths 1.25, 2.2, 3.5, 4.9, 12, 25, 60 
and 100 /im, and are labeled Bands 1-8 respectively. We obtained the 41 weekly- averaged DIRBE 
maps from the National Space Science Data Center (NSSDC), subtracted the zodiacal light model 
developed by the DIRBE team (Kelsall et al. 1998; Weiland et al, 1998, DIRBE Explanatory 
Supplement), and averaged the weeks together. 

The maps are stored in a cube format and pixelized into 6 faces of 256^ approximately 
square pixels ~ 0.3° on a side. For a finite beam the ideal map is convolved with the beam 
window function. For a circular top-hat beam similar to DIRBE's square top hat, the Fourier 
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transform of the window function is W{x) = [2Ji(x)/x]^, where Ji(x) is the first-order cyhndrical 
Bessel function. The measured power spectrum is then the product of the underlying power 
spectrum and the beam window function, i.e. -Pmeasured('7) = P2iQ)W{q'd) where i} reflects the 
beam size. On large scales (small q), W{x) ~ 1. Throughout the rest of the paper, P{q) will refer 
to the 2-dimensional power spectrum of the diffuse emission computed from DIRBE maps after 
deconvolution from the beam profile, i.e. P{q) = Pmcasuvedil) /W (q-d) . 

We determined the effective W and the effective beam size from maps of the beam response 
function archived at the NSSDC with 181x225 0.0065° square pixels, which we embedded in 256^, 
512^, and 1024^ pixel arrays. The results for Band 1 are shown in Fig. ^, and are similar for all 
eight bands. The solid line shows a top-hat beam profile with ■& = 0.4° which is slightly lower than 
the value of t9 = 0.45° adopted in Paper I. From eye fits we adopted = 0.4° in Band 1, 0.37° 
in Bands 2 and 3, and 0.35° in Bands 4 to 8. The combined effects of beam smearing from the 
motion of the beam during sampling, pixelization, and pointing error would increase the effective 
beam size by about 10%, and agree better with Paper I. The beam response function drops below 
10% at <0.15° and the pixelization prevents measurement below q~^ = 0.1°. 

4. Near-IR analysis 

4.1. Foregrounds and point source removal 

Foreground emission from the Galaxy and the solar system is the main problem in unveiling 
the expected CIB. At wavelengths less than 10 fim, the dominant foreground after removing the 
zodiacal light model is emission from stars in our Galaxy. Using the SAO Bright Star Catalog 
magnitudes and colors, we simulated J band maps. We found that aside from the large-scale 
shape of the galaxy and a few star clusters this foreground has an almost uncorrelated spatial 
distribution. At 12 and 25 /um (DIRBE Bands 5 and 6) the zodiacal dust is so bright that 
residuals from the model subtraction dominate the map structure outside the Galactic plane. At 
wavelengths of 60 and 100 fim (DIRBE bands 7 and 8) most foreground residual fluctuations come 
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from cirrus dust clouds in the Galaxy, and stars contribute little. 

The measured fluctuations contributed by point sources can be reduced substantially by 
identifying and removing them down to near the confusion noise limit. We used the point source 
finding routine developed by the DIRBE team and adopted in Paper I, and call it the "clipping" 
algorithm. The data from a selected region (patch) are first used to construct a smoothed model 
of the sky background emission for the whole patch, as follows. Surrounding each pixel, a window 
of size /size X /size is Searched for the minimal flux value. We use the minimal value rather that 
the median because the stellar brightness distribution is highly skewed. The map of this minimal 
(or lower envelope) flux is then fltted to a 2-dimensional polynomial surface of order ngt. This 
polynomial is in turn subtracted from the sky map, and the standard deviation cjo is calculated. 
Then, each of the pixels with flux above A'cutfo is masked out along with the 8 adjacent pixels, 
about twice the DIRBE beam area. This process of fltting a polynomial, identifying bright pixels, 
and masking them out, is iterated about 5 times, until no more fluxes above A'cut co found. In 
the near-IR no noticeable improvement is reached above ngt = 2. We discuss the dependence on 



/size in Sec. 4.5; for the results presented below we used /size = 5. The polynomial flt was used only 
for point source masking, no gradients were subtracted from the maps analyzed in this (near-IR) 
section, and the results at high galactic latitude (|fc| > 20°) were independent of this choice. 

For the near-IR, where foreground stars are the dominant source of fluctuations, there should 
be a clear correlation between the fluxes at different wavelengths. Furthermore, because most of 
the flux and fluctuations come from K and M stars, and the DIRBE pixels are large enough to 
include many stars, the dispersion in the color diagrams for the four near-IR bands should be 
small. In Paper I we used these properties to search for GIB fluctuations that have a different 
color due to redshifts or a different source spectrum. We construct a linear combination of maps 
at two different bands: 

Ai2 = SFi - f3uSF2. (17) 
Note that this A is unrelated to the A of Eq. ([l2[) . The variance of this map, Ca(0) = (Af2), is 
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minimized for 

- ~w¥r- ^^^^ 

Because the dispersion in the color index [5 = 5F1/6F2 is small, i.e. a/s = \/ — ^ 
with the average taken over all the pixels in the patch, the foreground contribution to C(0) will 
be reduced by a factor ~ (/3i2/<7/3)^ ^ 1. If the fluxes in adjacent bands do not correlate, the 
emission in the two bands comes from different sources or from measurement errors (noise). 

Since most of the predicted CIB comes from galaxies at 2; > 0.2 with typical stellar 
populations, its color should differ from that of the foreground stars, and the color subtracted 
maps should retain some of the CIB structure. If most of the CIB comes from high redshifts so 
that 

i^CIB = {SFi^ciB SF2,cib)/{SFIcib) < (19) 

the CIB fluctuations in the color subtracted maps will be larger than in the single band maps. We 
used Eq. (|l^) to make all-sky color-subtracted maps for all the adjacent band pairs, in order to 
search for a coordinate-independent part of the fluctuations. This method can be generalized to 
a multi-color subtraction method, e.g. minimize A123 = 5Fi — P6F2 — aSF^ with respect to a, (3. 
We applied this to the DIRBE maps, but without significantly different results. 



4.2. All sky variance analysis: C(0) 

We divided the sky into 384 patches, each with 32 x 32 pixels, and clipped each patch 
individually using the above procedure, with A^^^ut = 7, 5, 3.5 and 3. Because the foreground 
emission at these wavelengths is dominated by point sources, very few pixels are left for A'cut < 3. 
We performed the same analysis on 96 patches of 64 x 64 pixels, with similar results. 

Our star finding algorithm cuts deeper into the distribution than a simple interpretation 
might suggest. Each bright star is observed in about 5 pixels, depending on the position of the 
star. The noise distribution of these pixels in the absence of the bright star still has a non-Gaussian 
distribution with a non-zero mean. Therefore a star near the clipping threshold will be identified 
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in the pixel where the background noise fluctuations are greatest, rather than at the true star 
location. Hence, the clipping algorithm finds stars about 1 a fainter than A^^ut > as we confirmed 
with simulations. 

After clipping the 384 patches, we computed the color indices according to Eq. (18), the 
single band C(0), and the color-subtracted Ca(0) for each patch. Figure ^ shows histograms of 
the numbers of the remaining pixels after clipping. In the near-IR bands, about 350-450 pixels 
out of 1024, or about uncorrelated 75-90 beams per patch, remain for A'^cut = 3.5, 3, making the 
intrinsic uncertainty of C(0) less than 15%. 

In each patch there is a clear correlation between the fluxes in the adjacent bands. Paper I 
gave flux correlation plots for selected patches and an earlier model of the zodiacal light. The 
current data give similar results. For bands 1-4 all the patches have correlation coefficient R > 0.9 
between all the pairs of bands. The color indices have very small dispersion, af^/P < 10%. Because 
the dispersion in /? is so small, and most of the predicted CIB emission at these wavelengths comes 
from redshifts z > 0.1, the CIB fluctuations should not be removed by the color subtraction. 

Fig. 1^ plots the histograms of /? for A'cut = 3. The maps for the band pair [4, 5] have color 
index /? ~ because the dust in the solar system is not strongly correlated with the stars seen in 
Band 4. Between Bands 1 and 2, most of the patches have a color index of /3i2 — 2, typical of 
K-M giants. In band pairs [2, 3] and especially [3, 4] the spread in P is substantial, suggesting 
that star fluctuations are not the only source, or that their color is different in different regions. 
The range of color indices for the band differences [2 — 3] and [3 — 4] is reduced for \b\ > 20°, 
but nevertheless remains much wider than for the [1 — 2] maps. Most of the map peaks are due 
to stars, as we found them in simulated maps derived from the SAO Bright Star Catalog, using 
wavelength extrapolations appropriate to the tabulated spectral types. 

The lowest limits on C(0) for single bands are similar to those in Paper II, and for the color 
subtracted bands they are not very different from Paper I. However, we can now address the 
dependence of C(0) on Galactic coordinates. Fig. || shows a strong correlation of ^/CjO) with 
cosec|6|. Fig. |6| plots ^C(O) vs. cosec|6| for the color subtracted maps [1-2], [2-3], and [3-4]. 
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Although the color subtracted C(0) is a factor of ~ 10 lower than for single bands, the residual 
fluctuations still depend strongly on Galactic latitude. 



We would now like to extrapolate the data to estimate the extragalactic contribution to the 
fluctuations. We use analytic and numerical models to support a power-law dependence of the 
stellar contribution to C(0) on cosec|6|. We also show how to extrapolate the fluctuations due to 
Galactic stars to cosec|6| = 0, and that we expect no isotropic residual for C(0) due to stars. 

We assume that Galactic stars are distributed with a white noise power spectrum, i.e. they 
are spatially uncorrelated except for large scale distribution functions such as those in the Bahcall 
and Soneira (1980) and Wainscoat et al. (1992) models. Thus C(0) in the direction b is given by 
the following (Galactic) version of the Limber equation: 



Here the sum is taken over stars of type i with intrinsic luminosity Li, iv is the pixel solid angle, 
X = cosec|6|, r is the distance to each star, and rij is the number density of stars at that distance 
in the direction of b. The distance to the star is related to the height, Z, above the Galactic plane 
by r = Zcosec|6|. In our clipping algorithm, we remove sources exceeding the flux limit of A'cut 
X (standard deviation of the background polynomial fit). This standard deviation is approximately 
a from Eq. (^) above, if the integration limits are chosen to match the threshold for star detection 
and clipping. We find a very good correlation between the deviation of the fit and the value of the 
residual fluctuations in the patch; if all the histograms had the same non-Gaussian shape the two 
should be proportional. Hence, the lower limit radius can be approximated to depend on a itself 
via 



4.3. 



Near-IR fluctuations from Galactic stars 




(20) 



R^ = [Li/(47riVeutfT)] 



1/2 
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Equations (^), (^TJ) then allow us to determine a{b) for a given distribution of stars. 

We now derive the latitude dependence implied by these equations. As an example, we use 
a plane parallel exponential model for stars in the Galaxy: = no,iexp(— ^). This should be 
valid if our clipping radius R is small compared to the exponential scale length for the Galactic 
disk ~ (3 - 4) Kpc, and is in rough agreement with the disk in hydrostatic equilibrium with 



the Galactic gravity field (Mihalas & Binney 1981). Now Eq. (20) becomes 

'2 



'H^)-j:n,'^EA (22) 



where En{u) = t ^ exp{—ut)dt. We express (22) in terms of the total number of stars in the 
dRi interval: 

— = no^iUjRi exp(-— . (23) 
dRi hiX 

Then the equation for a{x) becomes: 

1 = ^cut E (^^^) [eMQ)E2{Q)]U=R,/h..- (24) 

Consider first the limit when Ri ^ hi. In this case g < 1 and the term in square 
brackets is near unity. The number of clipped stars per beam-width then becomes 
^> = ojJ2iJ(^"n{r;x)r'^dr ~ 3'^ J2i{RidNi / dR.;^ , leading to iV> = N~l/3. This is much 
smaller than we find in the data. Our clipping algorithm removes 9 pixels per star and cuts deeper 
than A'cut =3.5 would indicate, and leaves only about 40% of the pixels. This proves that q > 1, 
i.e. bright stars are seen out to the scale height of the disk and beyond. 

We now consider the variations in star counts with magnitude in a plane-parallel Galaxy. The 
differential counts in the magnitude interval dm are given by dlVfdTTi — co'[(ilnr^/(i77x]Tz(rj7^)r^, 
where oc IQ^-^"^ is the radial distance to the star of apparent magnitude m. In the limit of 

< /i, or n ~ constant, the counts converge to the uniform distribution limit, dN/dm oc lO*^'^™". 
In a plane-parallel Galaxy in which the radial structure can be neglected, the number density of 
stars depends only on \Z\ oc IqO-^™- gin \b\, and the differential counts in the direction x = cosec|6| 
can be related to those at the Galactic pole by 

dN ^ dN{m-5log^Qx) 

-J—\x = X Pole- (25j 

dm dm 
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Hence, once the counts at the pole are measured, we can evaluate the expected number counts in 
any direction b and then compute the fluctuations in the flux they produce via 



This formulation has two advantages: 1) the data on star counts in the polar regions can be used 
to evaluate the dependence of star fluctuations on cosec|6|; and 2) the only assumption is plane 
parallelism. We tested this assumption with both models and observations. 

We now determine dN/dm from the DIRBE data, and compare with prior observations. 
Figure ^ shows a Band 2 histogram of log^gdF — {F)\) for a 128^ pixel (38.4° x 38.4°) region at the 
North Galactic Pole (NGP). This logarithmic form is useful for displaying both large and small 
deviations, and we have normalized it to the equivalent number of stars per magnitude per square 
degree. The equivalent K magnitude is plotted on the top horizontal axis. With the DIRBE 
beam of 1.42 x 10~^ sr (4.5 pixels), F = 1 nWm^^sr^^ corresponds to rriK = 9.45. Because the 
GIB level and the smooth foreground contributions from dust and faint stars are not known in 
advance, we must subtract the mean of this confusion noise before seeking to measure brightnesses 
of individual bright stars. We define {F) as the mean flux of the map clipped with A'^cut = 3.5; in 
this case it is 66 nWm~^sr~^ or rriK = 5, near the turnover of the distribution. The total number 
of pixels with positive values of the flux after the subtraction is 12,359; the remaining 24.6% of the 
pixels have negative fluxes. We have plotted both positive and negative deviations to illustrate the 
asymmetry of the histogram far from the peak, as well as its symmetry near the peak. Negative 
deviations are statistical fluctuations, while very bright pixels are individual cataloged stars. For 
the clipped map at A^'cut = 3.5, a total of 6,380 pixels remain. Out of these, 2,982 pixels (47%) 
have positive fluxes, showing the symmetry of the flux distribution for the clipped map. 

Figure |^ also shows a Gaussian distribution with the variance cr ~ 10 nWm^^sr^^ computed 
with A'cut = 3.5, and for (7 = 9 nWm~^sr~^ . The difference between these two Gaussians is 
not large on this plot, showing the difficulty of direct detection of a sub-population of Gaussian 
fluctuations with dispersion of 5 nWm^^sr^^ . 




(26) 
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The NGP star counts were observed directly by Elias (1978). We show his data at = 1, 2.5, 
3.25 and 8 with N^/'^ error bars and our binning of his data. Further NGP data were obtained by 
the 2MASS survey in band, almost identical to the DIRBE Band 2, and were kindly provided 
to us by Tom Jarrett (1998). The cumulative counts from these measurements were shown in 
Fig. 1 of Beichman (1996) out to Kg > 15, who found that they follow dN /dm oc lO^ '^'" (cf. his 
Table 4). Actual 2MASS star counts from a region of 5 square degrees centered on the NGP are 
plotted in Fig. ^. The agreement between the DIRBE counts, the Elias (1978) and Jarrett (1998) 
data, and the B = 0.3 extrapolation is excellent over 15 magnitudes, or six decades in flux. South 
Galactic Pole counts from Minezaki et al. (1998) are also shown, and confirm the slope as well as 
the north-south symmetry of the Galaxy. At rriK < 1-5, the counts tend to the slope of B = 0.6 
coming from stars much closer than the scale height; if B were less than 0.6, the integrated star 
brightness would diverge at the bright end. 

The star counts agree with model predictions. Both Beichman (1996) in his Fig. 1 and 
Minezaki et al. (1998) in their Fig. 1 show that the counts are fitted well by extensions of either 
the Bahcall & Soneira (1980) or Wainscoat et al. (1992) models. An eyeball fit to their data gives 
B = 0.3 — 0.32 at K = 11. The Wainscoat et al. model at X = 11 shown in Fig. 1 of Minezaki et 
al. (1998) gives logdN/dm ~ 1.35, whereas continuation of the solid line in our Fig. |^ to ii' = 11 
gives logdN/dm = 1.3 if i? = 0.3 and 1.4 if B = 0.33. The agreement between the two slopes and 
normalizations is thus very good. Even the large-beam DIRBE instrument sees far beyond the 
scale height of the bright K band stars. 

We can now estimate the star brightness fluctuations. We begin with a power law for star 
counts in a plane-parallel Galaxy, dN/dm\po\e oc 10^"*. Then the fluctuation due to stars fainter 
than m at the Galactic pole becomes: 

a^lpoieoc r 10-°-S'"'+^"^'dm'. (27) 

J m 

We now derive the latitude scaling law implied by these assumptions. In our clipping method the 
lower limit on m in the integral is given by 10"'^-'^™ cx; A'cutc. Combining this with Eq. ( p5| ) we 



- 21 - 



find that the fluctuation in the total star flux at latitude b is given by the power law 

0- oc (cosec|6|)^-2. (28) 

A uniform, infinite star cloud has B = 0.6, so there is no latitude dependence in that case. If 
B = 0.4, then we expect a oc cosec|5|, while ii B = 0.3 then a oc cosec^|6|. The dependence of 



C(0) on cosec|6| can also be inverted via Eqs. (p5|), (26) to give a model for dN/dm. There is a 
unique relation between the two functions, for a plane-parallel Galaxy and our method of bright 
star extraction. 

The value of A'^cut enters indirectly, since it determines the magnitude m of the cutoff of 
detected stars and hence the magnitude at which B is important. Another measure of the star 
distribution function is the number of clipped stars, A>. In the power law case we have a very 
simple result: 

0.8-i? 2p+l 

- DA72 - OAT2 ' y^^) 



BNl, 3A2^, 



where p = {1.2/B) — 2 is the logarithmic slope in Eq. (p8[). The dependence of A> oc agrees 
with our results from DIRBE. Evaluating (^) for B = 0.6 (or < h) we find A> = l/3A^^j, in 
agreement with the exponential distribution of Eq. (24) with h = oo, while for B = 0.3 we find 



A^> = 5/3 A^y^-. The fraction of pixels clipped in the near-IR DIRBE maps at Acut = 3.5 is 



0. 



Because the clipping algorithm uses a mask with about twice the beam area, this corresponds to 
A> = 0.35 so the effective N^ut = 2.2 to 2.4 for p = 2.0 to 2.5. The effective A^ut is less than the 
input value of 3.5 because each star is observed in many pixels, each with its own noise. 



Similarly, Eqs. (|25|), (Pq) allow us to evaluate the expected variation of the residual fluctuation 

]^-|-0 S/B 1/3 

with Acut ) c" Of -^cut -III limit of a spatially uniform distribution, B = 0.6, and a oc A^y^ . 
For B = 0.3 it follows that a oc N^^^, and reducing Acut from 7 to 3.5 should decrease C(0) by a 
factor of 10, in agreement with the DIRBE data. 

To summarize this section, we have shown that, for dN/dm oc 10"^™ with B = 0.3 out to 
K < 15, and a plane-parallel Galaxy, we should recover a unique dependence for Galactic star 
fluctuations, ^C(O) oc (cosec|6|)^. This value of B is consistent with the Acut dependence of the 
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fraction of clipped pixels and the amplitude a of the remaining fluctuations in the DIRBE data. 



4.4. Extrapolating to cosec|&| = 

The scatter in Figs. |5| and ^ is large, because the plane-parallel model does not describe the 
large scale structure of the Galaxy. We divided the sample of 384 patches into four latitude bins 
and plot them against Galactic longitude I in Figs. |8| and ^. The dependence on longitude is 
significant for latitudes even as high as \b\ ~ 60° . However, for Galactic longitudes between 90° 
and 270° and high latitudes there is almost no longitude dependence. We therefore selected data 



with 90° < / < 270°, and plot the dependence of ^/C(0) on cosec|6| in Figs. 10 and 11. The scatter 
is much reduced. 

We tried several fitting functions for the 127 patches for which |6| > 20°, 90° < / < 270°: 



/Cfit(O) = a + yl(cosec|6|f , (30) 

where a presumably contains the cosmological and instrumental parts of the signal, 

Cfit(O) = a2 + A2(cosec|6|)2f, (31) 
, 2 

VCfit(O) = a + J2 Ai{cosec\b\y, (32) 



i=l 

and 



Cfit(O) = a exp(74cosec|6|), (33) 

minimizing, for example, ([(C(0) — Cfit)/C(0)]2) with respect to a,A,p. The fits for Eq. (31) are 
shown in Figs. 10 and ^ and in Table 1. Except for the color-subtracted [2 — 3] maps, all fits give 
positive values of a. The color-subtracted maps [2 — 3] and [3 — 4] have large scatter in the color 
indices at which the variance Ca(0) is minimized, and do not allow for a robust determination of 
a. Note that p > 1 as expected from the model. The star and GIB fiuctuations add in quadrature, 
as we demonstrated by simulation, so functional forms such as Eq. (^) are better justified than 
the other fits, although all give consistent results. The last row in Table 1 summarizes our limits 
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on a with errors corresponding to the extreme range from 92% confidence levels from all the fits, 



around our preferred central value from Eq. (31). 



We estimated the statistical uncertainty as follows. We define a relative variance 
o"o = min([(data — fit)/data]^) and a normalized XAr(a>^>p) = ([(data — fit) /data] ^)/(Tq. Fig. 
12 plots deviation histograms vs. ((5fit/c)^ where d^t = (C(0) — — A^a;^*')/-\/C(0). For purely 



Gaussian deviations the histograms would be straight lines of slope —1/2. 

We plot contours for Ax^ = 7 which for 3 parameters (a,p, A) corresponds to a confidence 



level of 92%. Fig. 13 shows thus determined confidence contours projected onto the (a,p) plane 
for Bands 1 to 4, according to Eq. (|3l|). The uncertainties shown in Table 1 correspond to the 
largest span of a in the panels in Fig. |l^; for any given value of "p the corresponding uncertainty 
levels are reduced significantly. This is the reason for the smaller formal uncertainties on a when 
an exponential fit - Eq. (|3^) - is assumed. Excluding from the analysis the patches that lie close to 



the Ecliptic plane further reduces the scatter in Figs. 10, O, but the contours are almost identical 



to those in Fig. 13. Likewise, keeping only the patches at 120° < / < 240° produces the same 
result. The results for iVcut = 7 and 5 are consistent but with a larger uncertainty, whereas for 
-^cut = 3 the results are close to those shown in Table 1. The values for p are in good agreement 
with the expectations from the star counts. For Band 2, the maximal spread is 1.5 < p < 2.6, 
which corresponds to 0.26 < B < 0.34 in good agreement with Fig. ^. 

To summarize, we have justified simple power law star count models and a plane parallel 
model of Galactic fiuctuations by comparing the predictions to the measurements. The model 
applies only far from the Galactic center and Galactic plane. There is a statistically significant 
residual term in the DIRBE data after extrapolation to zero cosec|6|, with consistent values for 
four different fitting functions and for different A'cut • Plate 1 of our Paper II shows visually that 
there is no obvious structure in the selected anti-center |6| > 20° data set. We checked the validity 
of our simple models by studying the residuals of the fits. Fig. 14 shows the residuals for Eq. (^) 



for each patch fiuctuations, ^/C(0) in Band 2. There is no apparent correlation with Galactic 
longitude or latitude and the residuals are similarly independent of the ecliptic coordinates. 
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To consider whether the instrument itself or the zodiacal light could be a source for this 
residual, we review the estimates in Paper I. There we constructed maps of the differences between 
different weeks of observation, and found noise levels ^/C{0) ~ 1.5, 0.3, 0.1 nWm^^sr^^ in Bands 
1 through 3 for A'cut = 3.5. These are much smaller than a in Table 1, and it is unlikely that a 
is due to noise. Similar arguments apply to errors in the zodiacal light modeling, which varies 
from week to week, and would contribute to the measurement error calculation. In any case, in 
the near-IR bands the contribution of the zodiacal light is small. It should be largest at 12 and 25 
fim, where the residual fluctuations are below 1 nWm^^sr^^ (Paper II and Sec. |^ of this paper), 
or < 1% of the total foreground. Therefore, the near-IR zodiacal modeling errors should be quite 
negligible compared to a. Furthermore, zodiacal light has a sharp cusp near the ecliptic plane, 
and zodiacal model errors would reflect this spatial dependence. We tried excluding patches in or 
near the ecliptic plane and found no change in a, p. The spatial correlation function that we find 



for A'^cut = 3 .5 (cf. Fig. 19) is significantly different from that expected for zodiacal light. The 
zodiacal light is smooth except near the ecliptic plane, where dust resonances and asteroid family 
debris are found. 

We conclude that the single band plots all indicate a positive and approximately isotropic 
residual term that is unlikely to be produced by either instrumental noise or errors in the 
zodiacal modeling. Since independent contributions add linearly to the combined variance, such 
a component would contribute only ~ 10% to the total dispersion of the confusion noise at the 
faint end of the flux distribution plotted in Fig. 0, and would not be detectable there. The Figure 
shows two Gaussians differing in dispersion by 10%, and they both seem reasonable fits. 



4.5. Analytic and Numerical Modeling 

Since the measured residual a can not be explained by known errors, we must investigate the 
stellar foreground fluctuations more carefully. To do so, we simulated the confusion noise process 
both analytically and numerically. 
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Our analytical approach assumes a power law for dN/dm, with a Fourier method to simulate 
the histogram of Fig. |^. Let V{F) be the probability distribution function to find a single star 
in the line of sight with flux F. The probability distribution function and its characteristic 
function G{f) are a pair of one dimensional Fourier transforms, G{f) = J V{F) exp(ifF)dF. The 
probability distribution function of finding two stars is the convolution of V with itself and the 
characteristic function for it is G'^{f). Similarly, for n stars the characteristic function is G^{f). 
For sufficiently large n the characteristic function tends to a Gaussian, which can be seen from 
expanding G{f) ~ 1 — /^(F^), (Chandrasekhar 1954). For a Poisson distribution with m stars per 
pixel on the average, the probability to find n stars in a given pixel is V{n) = exp(— m)/n!. 
Then the complex Fourier transform of the Poisson distribution of stars with an average of m 
stars per pixel is exp(m(G — 1)). For many stars and sufficiently small / or large fluctuations 
G{f) ~ 1 — ^/^(F^), and the previous expression converges to a Gaussian distribution whose 
width oc m~^/^. This prescription was implemented numerically assuming V{F) oc 10^^'^'^'°Sio ^ 
with B = 0.3. Convolution of the predicted star histogram with Gaussian measurement noise and 
Gaussian cosmic background fluctuations can be included by multiplying G{f) by a Gaussian. 
The effective beam size sets the maximum value of the apparent dN/dm in the confusion noise 
region. We are able to reproduce Fig. |^ very well, including the negative fluctuations and the 
transition from the confusion noise region to direct star detection. The result is robust in that the 
effective value of B measured at the bright end of the distribution is not altered by subtracting 
the mean value, and clipping with various A'cut has no effect on the value of B at the bright end. 
The apparent number of bright stars is about 10% larger than the input value, owing to the bias 
introduced by undetected stars of medium brightness near the detected bright stars. 

We also want to know whether there is any feature of our processing algorithm that could 
produce a spurious residual fluctuation a. For this, we need a numerical simulation of the sky, 
including both stars and possible cosmic terms. We used the DIRBE star model of Arendt et al. 
(1998), who implemented the 96-component star population model of Wainscoat et al. (1992), 
including spatial distribution models for the disk, bulge, halo, ring and arm populations. These 
were evaluated using the K-band dN/dm for the central pixel of each patch, and a 2-dimensional 
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uniform random number generator, and the stars out to X = 20 were placed at the centers of the 
pixels in the map. The flux from each star was then distributed among 9 pixels according to the 
measured or assumed beam shape, using delta-function, uniform, and Gaussian beam profiles. If 
the correct beam profile is used, the simulated maps processed with our algorithms match the 
DIRBE histograms in each patch very well. 

We also simulated CIB fluctuations, starting with Gaussian fluctuations in Fourier space with 
power spectrum P{k) = k^^'"^, as expected on the smallest DIRBE scales, and then multiplied by 
the top-hat beam window function determined in Sec. ^. The resultant field was then Fourier 
transformed to real space and normalized to the modeled cosmic variance, (Jsky We added these 
maps to the simulated star fluctuation maps and examined the results. 

As was discussed in Sec. |4.2| , our clipping algorithm clips effectively about ~ 1 a below 
the nominal value of A'cut ) because each star is seen in multiple pixels, each with its own noise 
fluctuations. Furthermore, according to Eq. (p9|), there should be no variation with cosec|6| in 
the number of stars clipped in each patch if the star counts follow the same B over the relevant 
range of magnitudes. To test this we computed the effective A''cut,eff = \J f ^> according to Eq. 
(|2^ ) assuming B = 0.3. The number of clipped stars, A^>, is half the number of clipped pixels 
since our mask has twice the beam area. The results are plotted against cosec|6| in Fig. |l5| for 
both the real sky and the model Galaxy. The data show that in the anti-center quadrants outside 
the Galactic plane, there is no trend with cosec|6|, and the effective values of A'^cut have very 
small dispersion. For the real sky, the mean in this range is A'^cut.cfT = 2.26 and the dispersion is 
cTcut = [(A'cut.eff) - (^cut,cff>^]^/^ = 0.06. For the model sky, the mean iV,ut,efr = 2.27, and the 
dispersion is 0.07. The effective Acut is independent of the simulated beam properties. 

We found a simple way to test the assumption of a plane parallel Galaxy. As Eq. (p5|) shows, 
the quantity x~'^{dN /dm)\x plotted against m — Slog^gX should be independent of x = cosec|6| 
and be equal to the counts at the Pole. We computed histograms for 14 patches of 64 x 64 pixels 
at \b\ > 20°, |/3cci| > 30° and 90° < / < 270°. The patches were clipped at A'cut = 3.5 and the 
average flux for the remaining pixels of each patch was subtracted from the map. The resultant 
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distribution of x^^dN/dm is plotted in Fig. for K band, with Poisson error bars for those points 
that contain at least 25 pixels (or 5.5 stars). The figure confirms that the stars are distributed in 
a plane-parallel way for this data set, with B ~ 0.3. The plot includes lines for the model Galaxy, 
with error bands. They agree with the observations within the errors except at the bright end, 
where the model is slightly low. 

We simulated skies with a cosmic term of a^-ky = 5 nWm^^sr^^ for the same 64 x 64 patches. 



The observations and the simulated model are plotted in Fig. 17, with fluxes measured in units of 



^C(O). If all the fluctuations are due to stars drawn from a power law distribution of the same 
slope, then all of the patches should follow the same line on this diagram. The DIRBE data and 
the synthetic Galaxy model look nearly identical, and there is no noticeable difference between 
Csky = 5 nWm^^sr^^ and 0. The DIRBE data and the models could be matched even better by 
adjusting the beam shape for the synthetic maps. The simulated confusion noise does not affect 
the amplitude and slope of the star counts at rriK < 5.5. 

Using our simulated models we also tested the scaling of the amplitude of the residual variance 
at given A'^cut as a function of the input isotropic component (Jsky This is important because our 
clipping algorithm might remove both cosmic and star fluctuations in a complicated and possibly 
non-linear way. We find that even at A'cut = 3.5, the residual variance C(0) is to good accuracy 
the sum of the variances from stars and the simulated GIB contribution. This shows that Eq. ( plj ) 
is better justified than the other choices for extrapolating to cosec|6| = 0. The choice of beam 
shape - Gaussian, delta-function, or uniform - has no effect on the effective A'^cut , but does lead to 
a systematic change in the residual C(0), with the largest beam area having the highest variance. 

Similarly, changing the size /size of the lower envelope region used by the point source 
recognition algorithm has no affect on A'^cut.cfr; but leads to a systematic dependence of 
C(0) oc l//sizc- Glearly, the larger the value of /size) the more diligently the algorithm recognizes 
point sources. We used /size = 5, corresponding to about 5 DIRBE beams. 

We can now predict what the Galactic star counts ought to be in order to reproduce the C(0) 
vs. X = cosec|6| relation found in the previous section. Gombining Eqs. (p5[) and (|2^) shows that 
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for a plane-parallel model the star contribution to the fluctuation is related to dNp/dm at the 
Pole via 

xa\x) = Fir lO-'-'y^dy, (34) 

where Fq is the zero magnitude flux; for Band 2 it is 6050 nWm~^sr~^ per DIRBE pixel. The 
lower magnitude for our clipping method is given by = — 2.51og]^o[/m-^cutc(x)/Fo]. Here fm 
is a factor accounting for the beam and the lower envelope used; based on the discussion in the 
previous paragraph we expect fm ~ 0.5 — 0.6. Differentiating both sides of Eq. ( |3^ ) leads to 



dNp ^ 0.4^ 0-4 In 10 d[xaHx)]/dx 

dm UNcutFo d[xM^)]/dx' ^ ' 

Here the right-hand-side should evaluated for x given by: 

m{x) = -2.51ogio[/mA^cut2:V(x)/Fo]. (36) 

Eqs. (|35D, (p^) form a closed set to determine the star counts required to reproduce the a{x). 
They thus provide an important consistency check between the measured effective A'cut > the 
C(0) - X relation, and the star counts at the Galactic Poles. They also show whether the 
cosec|6|-independent part of C(0) can be produced by the observed stars. 

We can measure the fm factor by normalizing the recovered star counts to the DIRBE data at 
rriK = 4, where the confusion noise is negligible, and find fm — 0.6. The solid line in Fig. |7| shows 



the recovered star counts at the Pole from C(0) given by Eq. (31) with a = with parameters 
{A^p) taken from Table 1. The agreement between the star counts inverted from the observed C(0) 
- X relation and the actual data is remarkably good, considering the statistical uncertainty in a 
and p, and given that all three (C(0) vs. x, A'^cut and dNp/dm) were determined independently. 
The dashed line in Fig. ^ shows the star counts required to reproduce Eq. (|30|) with a = 4.8 
nWm~^sr~^ according to Eqs. (^), (^). The line overshoots the data by several standard 
deviations at mpc > 5 and shows that this value of a cannot be produced by the observed Galactic 
stars. 

We can now use the simulated sky maps to test the extrapolations to cosec|6| = 0. We 
constructed simulated data for 384 patches of 32 x 32 pixels which contained both the K Band 
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Galaxy synthetic model, and a contribution from the CIB with o-gky varying from to 20 
nWm^^sr"-*^ . Processing the simulated maps with the standard algorithm, we find that the 
model without the CIB term is a good match to the DIRBE data at intermediate latitudes, but 
is significantly steeper as cosec|6| — > at |6| > 45°, and has a zero intercept within the statistical 
errors. On the other hand the simulated data curve at |6| > 45° flattens out for positive values of 

We also tried fitting the simulated sky fluctuations to the observations, using 
C(0) = ^mC^(csky) + CLm, and found the effective gain and offset Am,am- This approach 
has the advantage of including all that is known about the geometrical shape of the Galaxy, since 
the simulated maps use the Wainscoat et al. (1992) shapes for the disk, bulge, halo, ring and arm. 
If the model includes a cosmic term of the correct amplitude, then we should find A = 1, a„i = 0. 
We find that for a model agky = the value of = 1.4 it 0.3, consistent with unity, and the 
numbers for for positive values of crgky are similar. We could achieve ^ = 1 by better modeling 
of the DIRBE beam size, since we have already shown that relationship. 

We also confirm the values for the residual fluctuations a in Table 1. If the modeled sky has 
Csky = 0, we find y^a^ = 7.3l|;g nWm~^sr~^ in agreement with Table 1. Conversely, if we choose 
a model crgky similar to a in Table 1, then should be consistent with zero. We find — with 
the 92% confidence level for 2.8 nWm^^sr"-*^ < cTgky < 8.3 nWm^^sr"-*^ with the central value lying 
at cTgky = 5.2 nWm~^sr~^ . 

We can imagine only one possible feature in the DIRBE data processing that could lead 
to a false conclusion about the latitude dependence. We have shown that for a range of fixed 
beam size and shape, the sky models all predict a dependence of star fiuctuations on latitude 
that extrapolates to zero fiuctuations at zero cosec|6|. There is however some possibility that the 
effective beam size might depend on the line of sight and hence possibly on Galactic latitude. To 
explain the measurements, we would have to find an effect that systematically changes the effective 
beam area by an amount of the order of 30%, comparable to the fraction of the fiuctuations at the 
Galactic pole that seem to be of cosmic origin. 
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The DIRBE team measured the actual beam size very carefuhy from transits of bright stars 
through the beam, and did not find a significant dependence of beam profile on time or direction. 
However, as noted elsewhere, the effective beam profile for measuring fluctuations depends on 
several other effects. The effective beam area is approximately Q = '^'^{(^bcam/'^ + ^pixels + ^pointing)' 
where ^beam is the measured beam radius, ^pixels is the rms radius of the pixels, and ^pointing is 
the rms (vector) pointing error. The rms pixel size is dpix/S^^^ = 0.132° for a square pixel of 
side dpix, but for a rhombus of the same area and a 60° corner angle, it is increased by a factor of 
(2/^3)^/2 = 1.075 to 0.142°. I gnoring pointing error, this increases the effective beam area by 
3.2%, a negligible amount in this context. The measured DIRBE pointing error is 1.5 arcmin (1 
a), and increases the effective beam area by a fixed 3.9%. If, however, the pointing error were 
much larger than indicated by the statistics of the residuals from the pointing solution fits, and in 
addition were strongly dependent on Galactic latitude, the effect could be important for us. The 
pointing solution was the subject of extraordinary scrutiny by the COBE team, and such errors 
would have been noticed. We conclude that the latitude fitting method is not subject to errors 
due to changes of the effective beam size with latitude. 



4.6. Power spectrum 

Although the dominant spatial structure of the near IR maps is simply the white noise of 
stars, it is interesting to see whether large scale averages could reveal a GIB component. We 
computed both power spectra and angular correlation functions, and the results are shown in Figs. 



18 and 19. We describe the power spectra first. 



We used 96 patches of 64x64 pixels (or 19°xl9°), and computed spatial power spectra 
without star clipping, and with A'cut = 7 and 5, leaving over 90% and 80% of the pixels remaining 
at high latitudes. Smaller A^cut left too few pixels for reliable power spectra, and showed significant 
effects of the masking. We computed power spectra for single bands and for the color-subtracted 
maps, 6i — (352, with (3 evaluated for each patch according to Eq. (^). After each power spectrum 
was calculated, it was divided by the beam window function discussed in Sec. |^ to take out the 
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instrument signature. The fact that the power spectra are approximately flat confirms that stars 
are the dominant sources of fluctuations for A'^^ut = 5, and that the effective window function, 
accounting for pixeUzation, map distortion, beam smearing, and pointing errors, is correct. 

The single band power spectra also show large-scale gradients produced by the Galactic 
structure. The amplitude of the power spectrum decreases with A'cut) but the overall shape does 
not change appreciably, indicating that the beam mask effects are negligible for these values of 
A'cut • The amplitudes of the power spectrum for all angular scales are too large, and the shape 
too wrong, to allow for detection of the GIB structure. However, for some patches at A'cut = 5 the 
spectrum is close to that expected from Table 1. For A^^^ut of 5 or more, we can not extrapolate 
these power spectra to cosec|6| = for most of the scales probed by P{q), and for smaller A'cut 
there are too few pixels left. 



Fig. |18| shows the power spectrum for patch No. 7, at Galactic = (115°, 61°), and ecliptic 
(/?Ech Aeci) = (56°, 163°). In the single bands, the amplitude of the power spectrum decreases 
strongly after point source removal, but the shape remains approximately the same. The solid line 
in each panel shows the GIB signal according to Table 1, assuming that the GIB power spectrum 
has P oc g"^'^. This is a valid approximation for scales below one degree, but on large scales the 
power index of the GIB power spectrum may be different. The power spectrum of the foreground 
exceeds the estimated GIB by only a modest factor, particularly in Band 4. An instrument with a 
smaller beam might detect these fluctuations directly. 



Golor subtraction significantly reduces the foreground structure, as shown in Fig. 18. Without 
clipping, we find a reduction by a factor of ~ 30 between adjacent bands, although the spectra are 
still fiat from star fiuctuation noise. Galactic stars are dominant in Band 4 and zodiacal and cirrus 
emission are dominant in Band 5, so these bands are not strongly correlated. For this patch, color 
subtraction of Bands 4 and 5 shows structure dominated by cirrus dust, which has P (x with a 
steep n = —2 to —3; see Sec. The solid line in the lower panels shows the estimated small-scale 
GIB between Band 1 and Band 2 with logarithmic slope of —1.3 and C(0) taken from Table 1. In 
Bands [2-3], we plot y/C(&} =2 nWm-^gr-i . 
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For lower values of N^nt and q^^ > 1°, the amplitude of P{q) is comparable to that expected 
from the CIB according to Table 1. Furthermore, the slope of the power spectrum for this patch 
flattens significantly in the color-subtracted maps, and at some scales and bands approaches 
the logarithmic slope expected from the CIB. However, inspection of the power spectra shows 
that even in the color-subtracted maps, we do not detect a cosmological signal for A'^cut ^ 5. 
The spectra are all consistent with the white-noise stellar distribution and large-scale Galactic 
gradients. To check for a cosmological power spectrum consistent with Table 1 we would have to 
clip to much lower levels of N^ut where beam masking problems prevent reliable determination of 
the power spectrum. 

For A'^cut < 5, it is better to compute angular correlation functions, which are unbiased 
by masking effects, even for A'cut = 3.5 where less than 50% of pixels remain. We evaluated 
C{6) = {6F{x + 6)6F{x)), averaging over the remaining pairs of pixels separated by angular 
distance 9. (The points at different 6 are not statistically independent of each other. For A'cut = 
7 and 5 the correlation function is that of white noise, a delta function at zero lag, with a level 
consistent with the power spectrum plotted in Fig. It remains close to the zero-lag value on 
scales inside the beam, 9 < 0.5 and rapidly drops to very small (positive or negative) numbers on 
larger scales. For A'^cut = 3.5 for patches with low values of C(0) the correlation function flattens 
significantly, which would be consistent with CIB structure at the levels of Table 1. Fig. ^ shows 
the correlation function for Patch 7 at A'cut = 3.5, for 10 bins of 1 degree width. The straight line 
illustrates the slope of a CIB correlation function oc 6^^'"^ . The correlation function amplitude and 
slope are consistent with Table 1 and a CIB interpretation. Large-scale gradients in the star and 
dust populations are responsible for the positive correlation function at very large angles. 

A further comparison between the results in Table 1 and the power spectrum analysis can be 
made in the following way. For A'cut = 5 at each angular scale q~^, we evaluate the minimal value 
of P{q) for all the patches, and plot the value of a typical fluctuation, y/q'^P{q)/2TT, versus 6 = n/q 
in Fig. |2^. The slope of the data is roughly that of white noise, \/ q^P((i) j^i^ cx q^ implying that 
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we are still seeing Galactic stars. The upper limits over angular scales of 2° < 6 < 15° are: 

SFrmsie) < Ai-^)-\ (37) 

with A = 6, 2.5, 0.8, 0.5 nWm^^sr^^ for Bands 1 to 4 respectively. The process of finding an 
all-sky minimum of P2{q) produces a very smooth curve whose uncertainties are predominantly 
systematic. These are strong constraints on galaxy evolution. The shaded areas represent the 
power spectra according to Table 1, assuming that the CIB P{q) oc g"^ '^. Table 1 is consistent 
with the upper limits given by Eq. ([37|), except in Band 4. In Band 4, Table 1 would require a 
small scale for the turn-over in the CIB power spectrum, implying that much of it comes from 
high redshifts. 

5. Results in mid- to far-IR Bands 

5.1. Foregrounds and C(0) analysis 

At wavelengths greater than 10 fim, dust in the Solar system and the Galaxy produce most of 
the foreground emission. These sources are smooth on small scales, so do not necessarily prevent 
detection of CIB fluctuations. Odenwald, Newmark and Smoot (1998) detected < 100 nearby 
galaxies in the DIRBE data at wavelengths greater than 10 fim, but with the exception of M31 
and the Large and Small Magellanic clouds, the galaxies are unresolved. Because there are few 
discrete sources in the mid and far-IR DIRBE data, we can clip the DIRBE maps to lower values 
of A'cut and keep the same number of pixels as for the near-IR bands. 

Clipping to low values of A'cut would remove some fluctuations in the CIB as well as the 
foreground. If the CIB zero-lag signal were as high as the upper limit found in the maps, clipping 
at A'cut = 2, 2.5 or 3 levels would decrease the real C(0) by 20%, 10% and 0.001% respectively. 
Given that we find only upper limits, clipping down to A'cut = 2 is safe but lower levels would 
require interpretation. 

Fig. ^ shows histograms of the number of pixels remaining in each patch of 32x32 pixels 
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after removing point sources with N^ut = 2. Because of the extended nature of the foreground 
emission in these bands, we removed large-scale gradients first, using polynomials of order up to 
4. Fig. ^ shows the variation of the residual C(0) with cosec|6| at A^'cut = 2 for 90° < I < 270°. 
A significant fraction of the foreground emission still comes from the zodiacal light, even after 
subtraction of the DIRBE zodiacal light model, and contributes to the large scatter in the plots. 
Fig. also shows the dependence of C(0) on Galactic latitude for the 111 patches that also have 
Ecliptic latitude |/?eci| > 25°. The scatter is reduced but not enough to extrapolate to cosec|6|=0. 

In Bands 5 and 6 there is only a weak dependence on 6, and the amplitude of the typical 
fluctuations seen in these bands changes by only ~ 50% between cosec|6| = 1 and 3, showing 
approximate isotropy for \b\ > 20°. We interpret the lowest values of \/C(0) as upper limits on 
the CIB fluctuations. They are slight improvements over Paper II, and are shown in Table 2. 

In Bands 7 and 8 the dependence on Galactic latitude is more prominent. Extrapolation to 
cosec|6| = with Eq. (|30| ) gives values of a in agreement with the lowest ^/C{0) shown in Table 2, 
but with significant error bars. Therefore we again interpret the derived C(0) as upper limits on 
the CIB fluctuations. The slope of C(0) with Galactic latitude is consistent with a plane parallel 
Galaxy distribution, i.e. ^C(O) oc cosec|6|. The fluctuations are approximately proportional to 
the total brightness. 

The color diagrams at these wavelengths are not as clean as in the near-IR. In many patches 
no color correlations between the adjacent bands exist, and the color indices show large variations 
across the sky. Local variations in the parameters of the Galactic and solar system dust (e.g. 
density and temperature) are expected on large angular scales. No improvement in the CIB limits 
was achieved with mid-IR color subtraction. 

5.2. Mid- and Far-IR Power Spectra 

We computed power spectra for the same 96 patches of 64x64 pixels used for the near IR 
bands. Both single band power spectra and color-subtracted maps show a weak decrease with 
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clipping threshold, a consequence of the extended character of the foreground emission. The 
reduction of the foreground with the color subtraction is not as large as in the near-IR bands. In 
order to preserve information about P{q) at all angular scales, no gradients were removed before 
power spectrum analysis. 

The power spectra have a shape typical of the known cirrus distribution (cf. Gautier et al. 
1992). The spectrum P{q) is steep: P{q) oc g"" with n ~ (2.5 — 3). This is consistent with little 
small scale structure and, hence, no strong dependence in the resultant C(0) on A'^cut- It is also 
consistent with the power spectrum of cirrus emission measured on arc-minute and degree scales 



(Gautier et al. 1992; Wright 1998). Fig. |2^ shows power spectra for the same Patch 7 used in 
the near-IR analysis. It shows the steep power spectrum typical of the cirrus distribution with 
P{q) oc for most of the scales. Note that some of the fluctuations at 12 and 25 /im can be due 
to errors in modeling the zodiacal light. 

As in the near-IR analysis, we evaluated the minimal values of P{q) for each angular scale 
q~^ across the entire sky. Fig. ^ shows these minimal values of the fluctuation \/q'^P{q)/2TT as a 
function of angular scale ir/q. The minimal values come close to, but are not as small as, those 
in Table 2, which were evaluated at A'cut =2 after gradient subtraction. Even in the patches with 
the least fluctuations, the power spectrum is still as steep as P oc leaving the fluctuation 
y^q'-^P{q)/2Tr ~ const. This is indicative of cirrus emission in all areas of the sky. Infrared sky 
surveys with higher angular resolution should be able to reduce this contribution, and possibly 
uncover the GIB fluctuations. 



6. Conclusions 

Fig. ^ summarizes our results. At 1-5 /im, we find a positive residual fluctuation by 
extrapolating C(0) to zero cosec|6|. Based on detailed numerical and analytic models, this residual 
is not likely to originate from the Galaxy, our clipping algorithm, or instrumental noise. We 
conclude that this extra variance may result from structure in the CIB. The variance found in 



- 36 - 



this way from individual DIRBE Bands 1-4 is plotted with diamonds, with 92% uncertainty levels 



from Fig. 13. The results for color-subtracted maps are plotted with triangles at wavelengths 
halfway between the bands used. We find a positive residual in the color-subtracted map [1-2] 
between Bands 1 and 2, but not in other color subtracted maps. The color-subtracted map 
[1-2] has a unique color of /3i2 — 2 with little variation across the sky and the limit measures 
[{{SFi — 26F2)'^)]^^'^ . Taken at face value, these high values of the near-IR CIB fluctuations, if 
produced by evolving normal galaxy populations, would require substantial CIB fluxes. These 
would have to be above the estimates from the K-band galaxy counts, but are below the upper 
limits found by Hauser et al. (1998). The upper limits on CIB fluctuations at 10 - 100 /im are 
plotted with arrows, and are below 1 nWm~^sr~^ , with ^C(O) < 0.5 — 0.7 nWm~^sr~^ at 25 fim. 
These limits are lower than those in Paper II and imply strong constraints on how and when the 
early galaxies formed and evolved. 

On larger scales, 2° < < 15°, we obtain upper limits on the CIB fluctuations from the 
all-sky power-spectrum analysis: {9/5 ) x 6Fy-^s{9) < 6, 2.5, 0.8, 0.5 nWm^^sr^-'^ in Bands 1-4 
respectively. These limits, when taken in conjunction with our possible detection of the zero-lag 
CIB signal, limit the turn-over scale in the spectrum of the primordial density field to not much 
more than ~ 100/i~^Mpc. 

We fully recognize the difficulty of finding small fluctuations in the presence of larger 
fluctuations from foregrounds. While we have found no local explanation for our results, it is still 
quite possible that the fluctuations are not of cosmic origin, but come from some fault of the 
instrument, the data processing, or an unexpected feature of the Galactic foreground. The best 
way to resolve this uncertainty is to get better data, such as from a higher resolution sky survey 
with exceptionally good attention to flat flelding. This may be possible with satellites like SIRTF 
or the proposed NGSS, or even from 2MASS data or rocket data. 

We acknowledge support from the NASA Long Term Space Astrophysics grant 399-20-61-02. 
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and by Thomas Kelsall, for developing an exceptionally well calibrated and stable instrument, and 



-37- 



for producing public data archives with the complex zodiacal light foreground models removed. 
We have benefited greatly from conversations with Richard Arendt on applications of star models 
to the DIRBE data and for a careful reading and constructive comments on the manuscript of 
this paper. Our special thanks go to Tom Jarrett for providing the 2MASS star counts data. We 
thank the referee, Michael Vogeley, for careful reading of the paper. 



- 38 - 
REFERENCES 

Arendt, R. et al. 1998, Ap.J., 508,74. 

Bahcall, J. and Soneira, R., 1980, Ap.J.Suppl, 44, 73. 

Baugh, C. k Efstathiou, G. 1993, MNRAS, 265, 145. 

Beichman, C.A. 1996, in "Diffuse Infrared Radiation and the IRTS", eds. Okuda, H., Matsumoto, 
T. and Roellig, T., p.82. 

Beichman, C.A. & Helou, G. 1991, Ap.J., 370, LI. 

Boggess, N. et al. 1992, Ap.J., 397, 420. 

Bond, J.R., Carr, B.J. and Hogan, C. 1986, Ap.J., 306, 428. 

Chandrasekhar, S. 1954, in "Noise and stochastic processes", ed. N. Wax, Dover. 

Cowie, L. et al. 1994, Ap.J., 434 , 114. 

Djorgovski, G. et al. 1995, Ap.J., 438, L13. 

Dwek, E. et al. 1998, Ap.J., 508,106. 

Elias, J.H. 1978, A.J., 83, 791. 

Fall, S.M., Chariot, S. & Pei, Y.C. 1996, Ap.J., 464, L43. 

Fixsen, D.J., Dwek, E. , Mather, J.C., Bennett, C.L., & Shafer, R.A. 1998, Ap.J., 508, 123. 

Franceschini, A. et al. 1991, Ap.J.Suppl., 89, 285. 

Gautier, N. et al. 1992, A.J., 103, 1313. 

Hauser, M. et al. 1998, Ap.J. ,508, 25. 

Hacking, P.B. and Soifer, B.T. 1991, Ap.J., 367, L49. 

Jarrett, T. 1998, IPAC/Caltech, private communication. 

Jimenez, R. & Kashlinsky, A. 1999, Ap.J., 511, 16. ( |astro-ph/98023370 . 

Kashlinsky, A. 1998, Ap.J., 492, 1. 



-39- 

Kashlinsky, A., Mather, J., Odenwald, S. and Hauser, M. 1996, Ap.J., 470, 681. (Paper I) 

Kashlinsky, A., Mather, J., and Odenwald, S. 1996, Ap.J., 473, L9. (Paper II) 

Kelsall, T. et al. 1998, Ap.J., 508, 44. 

Loveday, S. et al. 1992, Ap.J., 390, 338. 

Lilly, S. et al. 1996, Ap.J., 460, LI. 

Lin, H. et al. 1996, Ap.J., 464, L60. 

Maddox, S. et al. 1990, MNRAS, 242, 43P. 

Malkan, M.A. and Stecker, F. 1998, Ap.J., 496, 13. 

Martin, C. and Bowyer, S. 1989, Ap.J., 338, 677. 

Mihalas, D. and Binney, J. 1981, "Galactic Astronomy", Freeman & Co., San Francisco. 
Minezaki, T. et al. 1998, A.J., 115, 229. 
Moore, B. et al. 1994, MNRAS, 269, 742. 

Odenwald, S., J. Newmark, & Smoot, G. 1998, Ap.J., 500, 554. 
Partridge, B. &; Peebles, P.J.E. 1967, Ap.J., 148, 377. 

Peebles, P.J.E. 1980, Large Scale Structure of the Universe, Princeton University Press. 

Picard, A. 1991, Ap.J., 368, L7. 

Saunders, W. et al. 1991, Nature, 349, 32. 

Saunders, W. et al. 1992, MNRAS, 258, 134. 

Shectman, S. 1973, Ap.J., 179, 681. 

Shectman, S. 1974, Ap.J., 188, 233. 

Shectman, S. et al. 1996, Ap.J., 470, 172. 

Silverberg, R.F. et al. 1993, Proc. SPIE, 2019, 180. 

Stecker, F., Puget, J-L., k Fazio, G.G. 1977, Ap.J., 214, L51. 



-40- 



Vogeley, M. 1998, ^stro-ph 97l"T20£ . 
Wainscoat, R.J., et al. 1992, Ap.J.Suppl., 83, 111. 
Wang, B. 1991, Ap.J., 374, 465. 
Weiland, J. et al. 1998, DIRBE Explanatory Supplement, |http:/ /www. gsfc.nasa.gov/astro/cobe/ 



dirbe_exsup.html 
Wright, E.L. 1998, Ap.J., 496, 1. 
Yoshii, Y. & Takahara, F. 1988, Ap.J., 326, 1. 



This manuscript was prepared with the A AS lATf^^X macros v4.0. 



-41 - 



Fit 


Band 1 
1.25 /Ltm 


Band 2 
2.2 iim. 


Band 3 
3.5 /Ltm 


Band 4 
4.9 /xm 


Band 1-2 


Band 2-3 


Band 3-4 


V'C(O) = a + AxP 
a 

P 


12.51^:? 

^ _-i_n 41 


2.181^:^1 


1 Q+0.6 

-■-^-0.7 

„ _i_n 45 

2-311^:45 


2.01°:^ 
2.79;S:°^ 


1.921^12 





1 4+^-2 

-■-•^-1.6 

1.47 


a 

P 


IS.SlylQ 


5.9l^j 


2 4+0.5 


9 n+0-25 
^•^-0.5 

-■O'^-O 78 


7 

'•D_2.4 

1 1 n+0 31 

-i--i'^-o:33 





1.23 


^C{0) = a + Y:t=iAix' 
a 


14.9±5.4 


5.6±2.2 


2.5±1.2 


2.7±1.2 


6.6 ±0.8 




0.9 ±0.9 


VC(0) = aexp{Ax) 
a 


9.7 ±0.6 


3.7 ±0.2 


1.3 ±0.1 


1.3 ±0.1 


5.3 ±0.2 


3.6 ±0.2 


1.5 ±0.1 


x= [\n^yC{0)-\na\/A 
a 


8.8 ±1.3 


3.3 ±0.3 


1.1 ±0.1 


0.8 ±0.1 


4.8 ± 0.6 


1.8 ±0.4 


1.2 ±0.1 


Summary 
















a 






o 4+0.5 


2.oi°:f 









Note: X = cosec|6| 



Table 1: Limits on a in nWm ^sr ^ ; the range corresponds to 92% confidence levels. 
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Band 5 
12 /xm 



1.0 



Band 6 
25 /xm 



0.5 



Band 7 
60 /xm 



0.8 



Band 8 
100 urn 



1.1 



Table 2: Upper limits on [C(0)]^/2 in nWrn-^sr"! for Bands 5-8. 
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Fig. |l| (a) Left: Linear scales of 0.5° in Eq. (|T2|) vs. z. The minimal scale is almost 
independent of cosmology, (b) Right: A{k) = (iij^^A;^i-*3(/c; 0))^''^ vs. k for APM spectrum. The 
CIB minimal relative fluctuations on the DIRBE beam scale are ~ A at the minimal value of the 
linear scale. 

Fig. ^: Window function and beam profile for DIRBE Band 1 beam. Plus signs show 
embedding in a 256^ pixel field, asterisks 512^ pixels, and diamonds 1024-^ pixels. Solid line is a 
top-hat profile with a beam radius of ■!? = 0.4°. 

Fig. ^: Histogram of pixels surviving clipping at A'cut = 3.5 in the near-IR DIRBE bands 
for 384 patches of 32 x 32 pixels. Solid line is Band 1, dotted Band 2, dashed Band 3, and 
dashed-dotted Band 4. 

Fig. ^: Histograms of color indices (3 from Eq. ( [I^ ) for 384 patches. Thin solid lines are 
all-sky; thick lines are patches with \b\ > 20°. 

Fig. |: y/C(fi} vs. cosec|6| for Bands 1-4 and iVcut = 3-5. 

Fig. ^: Same as Fig. |5| but for color subtracted maps. 

Fig. ^: K band clN/dm and DIRBE pixel histogram for 128^ pixels at NGP. Flux F is absolute 
value of deviation from mean of patch after clipping with A'cut = 3.5, measured in nWm~^sr~^ . 
Positive deviations are + signs, diamonds negatives. Poisson error bars for DIRBE data assume 
4.5 DIRBE pixels per star. Dotted lines are a Gaussian fit to data after clipping, and a 10% lower 
dispersion. Dash-dots are positive pixels remaining after clipping and dash-dot-dot-dot shows 
remaining negative pixels. Solid line shows counts inverted according to Eqs. (^5|) and (|3^) from 
Eq. ( |3l| ) without an isotropic component; long dashes are inverted counts if a comes from stars. 
Filled triangles are differential counts from Elias (1978) NGP measurements. Filled circles are 
differential NGP counts from 2MASS (Jarrett 1998). Open triangles are cumulative 2MASS counts 
multiplied by 0.3 In 10 to convert to differential counts for dN/dm oc lO'^'^™'. Filled diamonds with 
error bars are South Galactic Pole counts from Fig. 1 of Minezaki et al. (1998). 
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Fig. ^: Longitude dependence of ^/C(Wj for individual DIRBE bands at various Galactic 
latitudes. The increase of fluctuations towards the Galactic Center can be seen for any Galactic 
latitude. 

Fig. ^: Same as Fig. |8| only for color-subtracted maps. 



Fig. [ig: Plots of y/C(U) vs. cosec|6| for J, K, L, M bands for 90° < Z < 270°. Sohd lines are 
fits of Eq. using data for \b\ > 20°. 



Fig. 11: Same as Fig. Qd for color-subtracted maps. 



Fig. Histogram of fit residuals from Eq. (|0|) for Band 2 = [(C(0) - 
A2cosec2p|6|)/C(0)]2 in units of ctq. 



Fig. 13: 92% confidence limits on J, K, L, M and J-K for the fits to Eq. (31). The plus sign is 



the most likely value of {a,p) from Eq. (|3T|). 

Fig. 0: Scatter diagram for fit residuals (5fit = [VC(0) - y/C^JITj]/ y/C(&) for Eq. (|l]) for 
Band 2 versus / and b. 



Fig. 15: Effective clipping A'cut vs. cosec|6|, according to Eq. (|29|) for p = 2 or = 0.3. Left 



is DIRBE data and right is simulated Galaxy model. 



Fig. [1^: K-band DIRBE star counts in coordinates where a plane-parallel Galaxy would be 
a single line; x = cosec|6|. Data are for 64 x 64 pixel patches with \b\ > 20° and 90° < I < 270°, 
with Poisson errors shown for A^pix > 25. Confusion noise affects counts at rriK > 5.5. Lines show 
the model Galaxy: solid is mean x^'^dN/dm and dashes are the itl-sigma spread. 



Fig. 17: K-band star counts in the 64 x 64 patches outside the Galactic disk and away from 



the center, as a function of the absolute value of the flux deviation from the mean in units of 



y/C(0), for A'^cut = 3.5. Left shows DIRBE data and right shows the simulated Galaxy and a GIB 
fluctuation of dsky = 5 nWm~^sr~^ from Table ||. Plus signs are for negative F and diamonds for 
positive F; they overlap for small 



Fig. 18: Near IR power spectra for Patch 7, well above the Galactic and Ecliptic planes. Plus 



-45 - 



signs are before point source removal, asterisks show N^ut = 7, and diamonds show A'cut = 5. Sohd 
hnes show P{q) expected if the CIB power spectrum scales as q^^'^ for a in Table |l[ 



Fig. |19|: Angular correlation function evaluated for Patch 7 for A^cut = 3.5. The zero-lag value 
is plotted at ^ = 0.15°, and numbers below 0.1 nW^m~^sr~^ are not shown. The filled circle shows 
a from Table |l], normalizing a CIB correlation function with C{9) oc O^^'"^. Diamonds show the 
absolute value of the correlation function in the negative range. 



Fig. 2C: All-sky minimum fluctuation ^ q^P{cj)l2Tx versus vr/g. Shaded areas show the range 
expected from Table |l], assuming the CIB power spectrum scales as q~^'^ . 

Fig. Histogram of pixels surviving clipping at Acut = 2 in the mid- to far-IR DIRBE 
bands for 384 patches of 32 x 32 pixels. Solid line is Band 5, dots Band 6, dashes Band 7, and 
dash-dots Band 8. 



Fig. 0: vs. cosec|6| for Bands 5-8 and A'cut = 2. Top panels are 90° < / < 270°. 

Lower panels are 90° < Z < 270° and |/3eci| > 25°. 

Fig. p^: Upper panels are near-IR power spectra in for Patch 7. Lower panels are all-sky 



minimum fluctuations yJq^P{q) /2-k plotted vs. the scale T^/q. Plus signs are without clipping, 
asterisks show Acut = 7, and diamonds show A'cut = 5. 



Fig. 24: Summary. Diamonds are values for a with 92% uncertainties from Table 1. Triangles 
are for the color subtracted maps, shown at the mean wavelength for the two bands. The 1 — 2 
limits are the left triangle and 3 — 4 the right. Dashes with arrows are upper limits in Bands 4-8. 
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